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Abstract 

An  investigation  is  made  of  the  effects  of  the  change  of  the  porosity  of  the  gas  diffuser  layer  (GDL)  on  the  performance  of  a  proton  exchange 
membrane  (PEM)  fuel  cell.  The  analysis  of  fuel  cell  performance  with  non-uniform  porosity  of  GDL  is  a  necessity  because  the  presence 
of  liquid  water  in  the  GDL  leads  to  a  non-uniformly  distributed  porosity  in  the  GDL.  To  implement  this  performance  analysis,  a  half-cell 
model  which  considers  the  oxygen  mass  fraction  distribution  in  the  gas  channel,  the  GDL  and  the  catalyst  layer,  and  the  current  density  and 
the  membrane  phase  potential  in  the  catalyst  layer  and  the  membrane  is  investigated.  Four  continuous  functions  of  position  are  employed  to 
describe  the  porosity,  and  differential  equations  are  derived  based  on  oxygen  transportation  and  Ohm’s  law  for  proton  migration  and  solved 
numerically.  Results  show  that  a  fuel  cell  embedded  with  a  GDL  with  a  larger  averaged  porosity  will  consume  a  greater  amount  of  oxygen, 
so  that  a  higher  current  density  is  generated  and  a  better  fuel  cell  performance  is  obtained.  This  explains  partly  why  fuel  cell  performance 
deteriorates  significantly  as  the  cathode  is  flooded  with  water  (i.e.  to  give  a  lower  effective  porosity  in  the  GDL).  In  terms  of  the  system 
performance,  a  change  in  GDL  porosity  has  virtually  no  influence  on  the  level  of  polarization  when  the  current  density  is  medium  or  lower, 
but  exerts  a  significant  influence  when  the  current  density  is  high.  This  finding  supports  the  scenario  proposed  by  previous  studies  that  the 
polarization  at  high  current  density  corresponds  mainly  to  mass  transfer  through  (or  the  concentration  activation  of)  the  membrane  assembly. 
©  2002  Elsevier  Science  B.V.  All  rights  reserved. 
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1.  Introduction 

Proton  exchange  membrane  (PEM,  originally  referred  as 
solid  polymer  electrolyte)  fuel  cells  are  considered  to  be 
promising  and  environmentally  friendly  power  sources  for 
both  mobile  and  stationary  applications.  This  because  of  the 
attributes  of  high  energy  density  at  low  operating  tempera¬ 
ture  (70-90  °C),  good  performance  in  intermittent  operation, 
quick  start-up  capability,  zero  emissions,  and  minimal  prob¬ 
lems  from  component  corrosion  or  electrolyte  leakage.  Nev¬ 
ertheless,  polarization  occurs  during  operation  of  PEM  fuel 
cells  and  this  influences  dramatically  both  the  performance 
and  the  commercialization  of  the  technology.  Accordingly, 
the  development  of  a  theoretical  model  of  the  PEM  fuel  cell, 
as  well  as  corresponding  analyses,  are  crucial  to  gain  a  good 
understanding  of  the  effect  of  the  operating  conditions  on 
the  cell  potential  so  that  the  polarization  level  can  be  well 
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controlled.  An  extensive  review  of  the  development  of  such 
a  model  is  provided  in  [1-3]. 

Bernardi  [4]  studied  the  humidification  of  the  inlet  gas 
that  is  required  to  maintain  a  water  balance  in  a  PEM  fuel 
cell.  Verbrugge  and  Hill  [5,6]  developed  a  one-dimensional 
model  for  the  proton  exchange  membrane  (MEA)  based 
on  the  Nernst-Plank  equation.  Bernardi  and  Verbrugge  [7] 
employed  gas-phase  transport  and  membrane  models  to 
investigate  the  performance  of  a  gas-fed  porous  cathode 
bonded  to  a  polymer  electrolyte,  and  later  they  developed  a 
one-dimensional  isothermal  model  for  the  cell  [8].  Springer 
et  al.  [9,10]  also  derived  an  isothermal,  one-dimensional, 
steady- state  model  of  a  PEM  fuel  cell  to  simulate  perfor¬ 
mance  with  a  partially  dehydrated  membrane.  Following  the 
work  of  Bernardi  and  Verbrugge  [8],  Singh  et  al.  [11]  devel¬ 
oped  a  two-dimensional  isothermal  model  to  investigate  the 
effect  of  two-dimensionality  on  fuel  cell  performance.  Rowe 
and  Li  [12]  produced  a  one-dimensional  non-isothermal 
model  of  a  PEM  fuel  cell.  It  was  found  that  the  temperature 
distribution  within  the  cell  is  affected  by  the  phase  change 
of  water  in  the  electrodes,  and  that  the  partially  humidified 
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reactants  decrease  the  hydration  and  thus  increase  the  resis¬ 
tance  of  the  membrane.  To  consider  the  above-mentioned 
models,  however,  it  is  necessary  to  solve  the  coupled 
non-linear  differential  equations  by  complicated  numeri¬ 
cal  computations.  Also,  difficulties  may  be  encountered  in 
that  the  mass  transfer  limitations  in  the  gas  diffuser  layers 
(GDLs)  and  gas  channels  can  be  underestimated  [13]. 

To  avoid  these  difficulties,  a  simplifying  approach  has 
been  proposed  [13-15]  in  which  the  average  surface  over¬ 
potential  across  the  catalyst  layer  is  used  as  a  set  parame¬ 
ter  in  the  model.  This  approach  can  actually  decouple  the 
governing  equation  and  result  in  an  easier- to- solve  prob¬ 
lem,  so  that  a  more  complex  gas  channel  topology  or  a 
higher  dimension  of  the  cell  system  can  be  taken  into  ac¬ 
count.  A  good  example  can  be  seen  from  the  research  of 
Yi  and  Nguyen  [15],  who  used  this  approach  to  investigate 
multi-component  transport  in  the  porous  electrodes  of  PEM 
fuel  cells  using  the  interdigitated  gas  distributors.  The  model 
describes  the  two-dimensional  flow  patterns  and  the  distri¬ 
butions  of  the  gaseous  species  in  the  porous  electrode.  The 
model  predicts  the  current  density  generated  at  the  interface 
between  the  electrode  and  the  membrane  as  a  function  of 
various  operating  conditions  and  design  parameters.  It  was 
shown  that  the  averaged  current  density  generated  at  the 
positive  electrode  (cathode)  increases  with  higher  gas  flow 
rates,  with  a  thinner  electrode,  and  with  a  narrower  shoulder 
between  the  inlet  and  outlet  channels  of  the  interdigitated 
gas  distributor. 

Most  of  the  studies  mentioned  above  have  assumed  for 
simplicity  that  the  porosity  of  the  GDL  is  constant.  This, 
however,  may  not  reflect  the  importance  of  GDL  porosity 
on  fuel  cell  performance.  In  fact,  the  GDL  is  composed  of  a 
thin  layer  of  carbon  black  mixed  with  polytetrafluoroethy- 
lene  (PTLE),  a  carbon  paper  consisting  of  rigid  carbon  fib¬ 
bers  with  a  diameter  of  about  10  p,m,  and  a  porous  structure 
formed  by  carbon  fibbers  with  a  pore  size  of  20-50  jxm  and 
70-80%  porosity  of  volume.  Any  change  in  the  composi¬ 
tion  or  the  morphology  of  the  GDL  can  lead  to  a  substan¬ 
tial  influence  on  fuel  cell  performance  owing  to  a  porosity 
change  [16-19].  In  practice,  the  presence  of  water  generated 
in  the  GDL  during  fuel  cell  operation  can  also  change  the 
effective  porosity,  and  this  change  may  vary  with  time  and 
position.  Although  it  has  been  recognized  that  the  change 
of  content  and  morphology  of  the  GDL  can  influence  sig¬ 
nificantly  the  optimum  performance  of  the  fuel  cell,  rele¬ 
vant  studies,  especially  those  with  a  systematical  approach, 
are  still  rarely  reported  in  the  literature.  Accordingly,  to 
consider  the  effect  of  the  porosity  change  in  terms  of  the 
non-uniformity  of  porosity  in  GDL  on  fuel  cell  performance, 
it  is  necessary  for  work  to  be  done  for  a  precise  fuel  cell 
design. 

To  consider  the  non-uniformity  of  the  porosity  of  the 
GDL,  Gurau  et  al.  [20]  developed  a  one-dimensional 
half-cell  model  in  which  the  concept  of  the  so-called  effec¬ 
tive  porosity  was  employed  to  account  for  the  fact  that  the 
pores  may  be  partially  filled  with  liquid  water.  To  simplify 


the  problem  so  that  an  analytical  solution  could  be  obtained, 
the  authors  divided  the  GDL  into  a  series  of  parallel  layers, 
in  each  of  which  the  porosity  was  considered  to  be  uniform 
but  of  a  different  value.  The  results  showed  that  the  value  of 
the  limiting  current  density  decreased  when  the  GDL  was 
flooded  with  water  (a  smaller  porosity).  In  the  present  study, 
fuel  cell  performance  is  investigated  by  considering  the 
so-called  ‘half-cell  model’  [20];  special  attention  is  paid  to 
the  effect  of  non-uniform  porosity  on  fuel  cell  performance 
in  terms  of  physical  parameters  such  as  oxygen  consump¬ 
tion,  current  density,  power  density,  etc.  The  non-uniformity 
of  the  porosity  is  accounted  for  by  four  different  continuous 
functions  of  position,  each  of  which  has  a  different  averaged 
value  of  porosity  and  a  different  type  of  distribution  across 
the  GDL. 

In  the  following,  the  mathematical  model  and  corre¬ 
sponding  assumptions  are  given  in  Section  2,  in  which  the 
governing  equations  and  associated  boundary  conditions 
for  both  the  oxygen  mass  fraction  and  the  membrane  phase 
potential  are  presented.  Mathematical  models  of  porosity 
are  also  discussed.  In  Section  3,  the  computational  re¬ 
sults  are  presented  to  illustrate  the  effect  of  non-uniform 
porosity  on  the  distribution  of  the  oxygen  mass  fraction 
across  the  GDL  and  the  catalyst,  as  well  as  on  the  current 
density  and  membrane  phase  potential  across  the  catalyst 
layer  and  the  membrane.  In  Section  4,  the  performance 
of  the  fuel  cell  in  terms  of  the  polarization  curve  is  dis¬ 
cussed  by  again  focusing  on  the  effect  of  porosity.  It  is 
also  shown  that  the  present  numerical  results  serve  well 
to  explain  the  physical  meaning  of  the  performance  model 
proposed  previously.  Linally,  in  Section  5,  several  con¬ 
cluding  remarks  are  made,  where  the  relation  between  fuel 
cell  performance  and  the  porosity  of  the  GDL  is  explicitly 
identified. 


2.  Mathematical  model  and  assumptions 

A  schematic  of  the  positive  electrode  (cathode)  side  of  a 
PEM  fuel  cell  is  shown  in  Lig.  1.  This  consists  of  four  dif¬ 
ferent  regions,  namely,  the  gas  channel,  the  GDL,  the  cata¬ 
lyst  layer,  and  the  membrane.  In  the  GDL  (usually  a  carbon 
fiber  paper),  the  morphology  varies  widely,  as  determined 
by  the  proportions  of  the  constituents.  Moreover,  in  order  to 
consider  the  presence  of  the  water  in  the  GDL,  the  effective 
porosity  is  taken  into  account.  Since  the  presence  of  liquid 
water  is  intrinsically  non-uniform  in  the  GDL,  the  porosity 
of  the  GDL  will  be  a  function  of  position,  instead  of  a  uni¬ 
form  porosity  or  a  piece-wise  uniform  model  as  considered 
by  previous  studies,  see  for  example  [14,20].  In  the  present 
paper,  four  different  models  of  the  porosity  of  the  GDL  are 
considered,  and  are  described  by  four  different  mathemati¬ 
cally  continuous  functions  (Lig.  2). 

(1)  Model  1  represents  constant  porosity  and  the  mean 
porosity  of  the  GDL  (x)  =  0.4  is  considered. 
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(2)  Model  2  is  a  linear  function  described  by  £d(x)  = 
—  1.143  x  10_3v  +  0.7,  which  is  a  mean  porosity  equal 
to  0.523. 

(3)  Model  3  is  a  convex  exponential  function  described  by 
ed(x)  =  exp  [-7.5  x  I0~6x2  -  2.04  x  I0~4x  -  0.357J, 
which  is  a  mean  porosity  equal  to  0.551. 

(4)  Model  4  is  a  concave  exponential  function  described  by 
ed(x)  =  exp  [7.5  x  10“6x2  -  5.05  x  10“2jc  -  0.357J, 
which  is  a  mean  porosity  equal  to  0.472. 


0  70  140  210  280  350 
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Fig.  2.  Porosity  distribution  function  of  four  models  considered. 


Model  1  of  uniform  porosity  allows  a  comparison  to  be 
made  between  the  present  numerical  results  and  the  ana¬ 
lytical  results  of  Gurau  et  al.  [20].  The  other  three  mod¬ 
els  are  chosen  to  account  for  the  effective  porosity  influ¬ 
enced  by  the  presence  of  water  in  the  GDL  during  fuel  cell 
operation.  The  use  of  an  exponential  function  is  primar¬ 
ily  due  to  the  convenience  of  the  mathematical  derivation. 
There  are  no  existing  experimental  data  by  which  the  real 
change  of  porosity  with  the  presence  of  water  in  the  GDL 
can  be  confirmed.  Accordingly,  no  model  is  thought  to  be 
particularly  suitable  for  a  cell  under  a  certain  operational 
condition.  The  present  exponential  models  are  therefore  be¬ 
lieved  to  be  able  to  account  for  real  conditions  to  a  certain 
extent. 

The  mathematical  model  of  the  present  system  is  made 
on  the  basis  of  several  assumptions  [20]:  (a)  the  model  is 
one-dimensional;  (b)  the  system  is  isothermal  and  steady;  (c) 
the  air  is  taken  as  an  ideal  gas;  (d)  the  air  is  fully  saturated 
with  water  vapor;  (e)  the  membrane  is  maintained  in  a  fully 
hydrated  condition;  (f)  oxygen  is  transported  to  the  catalyst 
site  as  a  gaseous  component;  (g)  diffusion  is  the  only  mode 
of  transport. 

2.7.  Governing  equations  of  oxygen  mass  fraction 

The  species  are  transported  through  the  gas  channel,  the 
GDL  and  the  catalyst  layer  only  by  diffusion.  The  governing 
equation  for  the  oxygen  mass  fraction  in  the  gas  channel 
is: 

d  (  dTg\ 

dH^df)=0’  (1) 

where  p  is  the  density  of  mixed  gases,  Do2  the  diffusion 
coefficient  of  oxygen,  Tg  the  oxygen  mass  fraction  in  the 
gas  channel  and  v  the  position.  In  the  GDL,  the  equation  is 

^  {[edWfVfo^J  =0,  Lg<x<Ld  (2) 

where  £d(x)  and  ra  are  the  porosity  function  and  tortuosity 
factor  of  the  GDL,  respectively,  and  Y&  is  the  oxygen  mass 
fraction  in  the  GDL.  Since  the  GDL  is  a  porous  medium,  a 
Bruggeman-type  correction  [21]  can  be  applied  to  account 
for  the  effective  diffusion  coefficient.  In  the  catalyst  layer, 
the  oxygen  is  consumed  and  the  governing  equation  is: 

dx  \  c  PDo2~dx  )  =  4FM°2  L**X*L'  (3) 

where  sc  and  tc  are  the  uniform  porosity  and  tortuosity 
factor  of  the  catalyst  layer,  respectively,  Yc  is  the  oxygen 
mass  fraction  in  the  catalyst  layer,  F  =  96,487  the  Fara¬ 
day  constant,  Mo2  =  32  the  molecular  mass  of  oxygen, 
j  =  2 FkYc  exp[2a Fq/RT]  the  cathode  transfer  current  den¬ 
sity,  in  which  k  is  the  reaction  rate  constant,  a  the  transfer 
coefficient,  r\  the  surface  overpotential,  R  =  8.314  is  the 
universal  gas  constant,  and  T  is  the  absolute  temperature. 
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The  boundary  conditions  at  different  positions  accounted 
for  by  the  oxygen  mass  fraction  are: 

7g  =  F0 

dFg 

lg  =  *d. 


,  v  =  0 

(4a) 

r  did 

=  [£<i(x)]T*Do2-f,  x  =  L% 

(4b) 

,  x  =  Lg 

(4c) 

did  t  dFc 

Do2  -j—  =  e?  do2  ~r~ ,  x  =  Ld 
dx  dx 

(4d) 

x  =  Ld 

(4e) 

x  =  Lc 

(4f) 

d Yc 
dx 

The  density  p  is  calculated  from  the  ideal  gas  law  [21],  i.e. 
PM 


P  = 


RT 


(5) 


where  P  is  the  pressure  and  M  the  molecular  mass.  The 
molecular  mass  of  the  gas  mixture  is  obtained  by: 


M  = 


1 


(Yi/MoJ  +  (Tw/Mw)  +  ((1  -  Yw  -  Yi)/MN2 )  ’ 


i  =  g,  d,  c 


(6) 


where  Yw  is  the  water  vapor  mass  fraction,  the  subscripts 
O2,  w  and  N2  mean  oxygen,  water  vapor  and  nitrogen,  re¬ 
spectively.  The  mass  fraction  of  the  saturated  water  vapor  is 
given  by: 


1  w  — 


PM 


(7) 


where  the  water  saturation  pressure  P^at  is  calculated  from 
the  following  expression  [22]: 

In  .(P“‘)  =  -  5800^2206  +  1.3914993 

—  0.048640239(7’ +0.41764768  x  10“4)P2 

-  0. 14452093  x  10_7(P3+6.5459673)ln(P)  (8) 

The  oxygen  diffusion  coefficient  is  a  function  of  the  binary 
diffusion  coefficients,  namely: 


1  ~Yi 


1 


°°2  (Fw/(Do2-wMw))+((l-Fw-y,)/(Do2-N2MN2))  M’ 

i  =  g,  d,  c  (9) 

The  binary  diffusion  coefficients  Do2-w  and  Do2-n2  are  then 
calculated  by  using  the  Slattery  and  Bird  estimation  [21]: 


2.2.  Governing  equations  of  membrane  phase 
potential 

The  membrane  phase  potential  exists  in  the  membrane 
and  the  catalyst  layer,  and  accounts  for  the  difference  be¬ 
tween  the  surface  overpotential  and  the  potential  of  electron¬ 
ically  conductive  carbon.  In  the  present  system,  the  averaged 
surface  overpotential  across  the  catalyst  layer  is  considered 
constant,  and  is  used  as  a  parameter  of  the  model.  In  the 
membrane,  the  governing  equation  of  the  phase  potential  is 
given  as: 


d 

dx 


d0m 

r - 

dx 


=  0,  Lc  <  x  <  Ln 


(ID 


where  o  is  the  ionic  conductivity  of  the  membrane  and  <Pm 
the  membrane  phase  potential  in  membrane.  In  the  catalyst 
layer,  the  equation  becomes: 


d 

dx 


( v  dM 

£  loncr -  = 

V  d  x  ) 


,  I  —  j >  —  %  —  Lq 


(12) 


where  £  and  T{on  are  the  porosity  and  tortuosity  of  the 
ionomer  in  the  catalyst  layer,  respectively,  <PC  the  membrane 
phase  potential  in  the  catalyst  layer. 

The  boundary  conditions  applied  at  different  positions 


are: 

^=0, 

dx 

x  =  Ld 

(13a) 

g^ion  d^C  _ 

d^m 

(13b) 

dx  “ 

,  X  —  Lc 
dx 

o 

II 

O 

3 

x  =  Lc 

(13c) 

0, 


(13d) 


In  above  equations,  the  membrane  ionic  conductivity  a  can 
be  obtained  by  the  empirical  expression  [9] : 

1  1 
, 303  ~  7\ 


0T=  (0.5139A.  —  0.326)exp 


1268 


(14) 


where  X  is  the  membrane  water  content  at  the  interface.  A 
typical  function  form  of  X  measured  for  Nation  117  can  be 
found  in  [9],  and  has  a  value  between  14  and  22,  as  dictated 
by  the  temperature  of  water. 

2.3.  Governing  equations  of  current  density 

In  the  membrane  and  the  catalyst  layer,  the  current  den¬ 
sity  is  governed  by  Ohm’s  law  for  proton  migration.  The 


J\,m  = 


10.13  x  fl(7yyrcr,k7’cr,m)fc(Pcr,kPcr,m)1/3(7’cr,k7’cr,m)5/12((l/Mk)  +  (1  /Mm))1'2 


(10) 


where  a ,  b  are  the  constants  in  the  Slattery  and  Bird  estima¬ 
tion  for  the  binary  diffusion  coefficients  and  the  subscript 
cr  accounts  for  the  critical  values. 


equations  are  written  as: 

r.  d0c 

ic  =  £  loncr— — ,  Ld  <x  <  Lc 

dx 


(15) 
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where  i  is  the  current  density;  the  subscripts  c  and  m  mean 
catalyst  layer  and  membrane,  respectively,  After  obtaining 
the  values  of  the  membrane  phase  potential,  the  calculation 
for  current  density  can  be  carried  out  easily  from  Eqs.  (15) 
and  (16)  without  employing  boundary  conditions.  The  above 
equations  are  the  same  as  those  used  by  Gurau  et  al.  [20], 
except  that  the  porosity  of  the  GDL  is  considered  as  a  func¬ 
tion  of  position,  as  shown  in  the  beginning  of  Section  2. 

3.  Effects  of  change  of  porosity 

The  resulting  steady-state,  non-linear,  differential  equa¬ 
tions  shown  in  the  previous  section  are  solved  by  a  fi¬ 
nite  difference  numerical  approach.  The  central  differencing 
scheme  is  applied,  which  is  a  routine  procedure  and  is  de¬ 
scribed  many  numerical  textbooks.  Thus,  it  is  not  necessary 
to  show  explicitly  the  discretized  equations.  The  base  case 
conditions  are  the  same  as  those  used  by  Gurau  et  al.  [20]; 
the  reader  is  referred  to  their  Table  1  for  the  details. 

3.1.  Effects  on  oxygen  mass  fraction 

The  effects  of  porosity  change  on  the  oxygen  mass  frac¬ 
tion  across  the  gas  channel  and  the  GDL  catalyst  layer  are 
examined  for  three  different  surface  overpotentials,  namely, 
y)  —  0.25,  0.28  and  0.36.  The  computational  results  are 
given  in  Fig.  3.  Results  show  that,  for  the  same  porosity  in 
the  GDL,  a  higher  surface  overpotential  r\  leads  to  a  lower 
oxygen  mass  fraction  in  the  gas  channel.  This  is  due  to  the 


Fig.  3.  Effects  of  porosity  change  on  distribution  of  oxygen  mass  fraction 
across  gas  channel,  GDL  and  catalyst  layer  for  three  different  surface 
overpotentials. 


fact  that  a  larger  surface  overpotential,  77,  corresponds  to  a 
more  active  electrochemical  reaction,  which  leads  to  an  in¬ 
crease  in  the  oxygen  consumption  in  the  gas  channel.  On 
the  other  hand,  for  the  same  overpotential,  a  higher  aver¬ 
aged  porosity  of  the  GDL  results  in  a  smaller  oxygen  mass 
fraction  in  the  gas  channel.  This  is  because  a  larger  poros¬ 
ity  of  the  GDL  allows  a  higher  consumption  of  oxygen  in 
the  catalyst  layer.  For  the  same  reason,  when  the  overpo¬ 
tential  is  smaller,  the  consumption  of  oxygen  is  smaller  and 
the  oxygen  consumption  is  significantly  less  influenced  by 
the  porosity  change  in  the  GDL.  See,  for  example,  the  case 
of  77  —  0.25,  in  which  the  oxygen  mass  fractions  in  the 
gas  channel  given  by  the  four  different  porosity  models  are 
virtually  the  same.  Under  such  circumstances,  the  electro¬ 
chemical  reaction  is  dominated  by  the  surface  kinetics  in  the 
catalyst  layer. 

In  the  GDL,  the  effect  of  a  continuous  change  of  poros¬ 
ity  can  be  seen.  When  r\  —  0.36,  for  example,  the  oxygen 
mass  fraction  changes  with  position  and  follows  closely  the 
porosity  distribution  across  the  GDL,  i.e.  as  stated  previ¬ 
ously,  the  oxygen  mass  fraction  is  smaller  when  the  poros¬ 
ity  is  larger,  and  vice  versa.  To  some  extent,  this  reflects 
the  fact  that  the  oxygen  mass  fraction  is  not  uniform  across 
the  GDL,  especially  in  the  region  adjacent  to  the  inter¬ 
face  between  the  GDL  and  the  catalyst  layer.  It  is  inter¬ 
esting  to  note  that  the  difference  in  oxygen  mass  fraction 
due  to  the  porosity  change  in  the  GDL  decreases  as  77  be¬ 
comes  larger  because,  for  a  larger  77,  the  electrochemical  re¬ 
action  is  more  active  and  more  oxygen  is  consumed.  When 
77  —  0.36,  the  electrochemical  reaction  is  so  strong  that  the 
oxygen  is  almost  depleted  completely  in  the  catalyst  layer, 
and  differences  in  the  oxygen  mass  fractions  between  the 
four  models  disappear.  For  low  77,  the  difference  remains 
because  under  such  circumstances,  the  electrochemical  re¬ 
action  is  dominated  by  the  surface  kinetics,  so  that  less 
oxygen  in  the  catalyst  layer  is  consumed  and  the  amount 
of  oxygen  remaining  increases  as  the  porosity  decreases. 
Consequently,  the  difference  becomes  more  significant  as  77 
decreases. 

The  change  in  the  oxygen  mass  fraction  across  the  cata¬ 
lyst  layer  due  to  the  porosity  change  in  the  GDL  for  77  = 
0.25,  0.28,  0.36  is  illustrated  in  Fig.  4.  Note  that,  because  the 
catalyst  layer  is  very  thin  and  the  electrochemical  reaction 
within  the  layer  prevails,  the  variation  in  the  oxygen  mass 
fraction  across  the  layer  is  not  at  all  significant  for  all  three 
values  of  77,  and  cannot  be  detected  on  the  decimal  scale. 
On  a  logarithmic  scale,  however,  the  variation  in  the  case  of 
77  =  0.36  can  be  clearly  seen,  i.e.  the  oxygen  mass  fraction 
decreases  exponentially  from  the  GDL  to  the  membrane.  Al¬ 
though  all  the  variations  are  very  small,  the  averaged  oxy¬ 
gen  mass  fraction  of  77  =  0.36  is  the  smallest  because  the 
oxygen  mass  fraction  decreases  as  77  increases  due  to  the 
dominance  of  the  electrochemical  reaction.  By  contrast,  the 
effect  of  the  porosity  change  in  the  GDL  on  the  oxygen  mass 
fraction  is  insignificant  for  all  the  three  values  of  77.  This  is 
because  virtually  all  the  oxygen  has  been  consumed  within 
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Fig.  4.  Variations  of  oxygen  mass  fraction  in  catalyst  layer  for  three 
different  surface  overpotentials  and  four  different  porosity  models. 

the  catalyst  layer.  Thus,  the  differences  between  the  oxygen 
mass  fractions  in  the  four  porosity  models  are  insignificant. 

3.2.  Effects  on  current  density 

The  current  density  distributions  in  the  catalyst-layer- 
membrane  domain  of  the  four  porosity  models  of  the  GDL 
when  r)  —  0.25,  0.28,  0.36  are  given  in  Fig.  5.  The  current 
density  increases  significantly  in  the  catalyst  layer  and  ap¬ 
proaches  to  the  limiting  value  in  the  membrane,  which  in¬ 


Fig.  5.  Effects  of  porosity  change  on  distribution  of  current  density  across 
catalyst  layer  and  membrane  for  three  different  surface  overpotentials. 


creases  with  r\  because  a  higher  overpotential  results  in  a 
more  active  electrochemical  reaction  and  thus  a  higher  resul¬ 
tant  electronic  current.  For  the  same  r)  but  different  porosity 
model  of  the  GDL,  the  current  density  increases  with  the  av¬ 
eraged  porosity.  Again,  this  is  obviously  due  to  the  fact  that 
a  larger  porosity  in  the  GDL  leads  to  a  larger  consumption  of 
oxygen  in  the  catalyst  layer,  and  thus  a  larger  current  density 
is  generated.  The  difference  between  the  highest  and  lowest 
current  density  can  be  as  much  as  16%  when  r)  —  0.36  and 
the  averaged  porosity  differs  by  about  30%  (the  difference 
between  models  3  and  4).  The  effect  of  variable  porosity  on 
the  current  density  is  not  seen  at  all  in  Fig.  5.  This  is  because 
the  electronic  current  occurs  only  in  the  catalyst  layer,  in 
which  the  chemical  reaction  is  influenced  only  by  the  amount 
of  oxygen  transported  from  the  GDL  into  the  catalyst  layer, 
and  is  independent  of  the  porosity  distribution  in  the  GDL. 

3.3.  Effects  on  membrane  phase  potential 

The  variation  of  the  membrane  phase  potential  0  in  the 
catalyst-layer-membrane  domain  for  the  different  porosity 
models  and  surface  overpotential  r\  =  0.25,  028,  0.36  is  il¬ 
lustrated  in  Fig.  6.  The  membrane  phase  potential  is  lower 
when  the  porosity  of  the  GDL  is  higher.  For  instance,  model 

3  is  for  the  highest  averaged  porosity  but  the  lowest  phase 
potential  among  the  four  models  considered;  whereas  model 

4  is  for  the  lowest  averaged  porosity  but  highest  phase  po¬ 
tential.  The  difference  between  the  highest  and  lowest  mem¬ 
brane  phase  potential  is  24%  when  r)  —  0.36.  A  higher 
porosity  in  the  GDL  leads  to  a  larger  current  density,  and 
thus  a  higher  loss  of  the  membrane  phase  potential.  As  a 
result,  model  3  (i.e.  highest  averaged  porosity,  lowest  mem¬ 
brane  phase  potential)  represents  the  best  cell  performance 


Fig.  6.  Effects  of  porosity  change  on  membrane  phase  potential  across 
catalyst  layer  and  membrane  for  three  different  surface  overpotentials. 
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among  the  four  porosity  models.  For  the  same  physical  rea¬ 
sons,  when  the  same  model  of  porosity  is  considered,  the 
loss  of  membrane  phase  potential  is  larger  when  the  overpo¬ 
tential  is  larger,  which  again  results  in  a  better  performance 
of  the  fuel  cell. 


4.  Physical  explanation  of  performance  model 

Models  to  describe  the  steady-state  performance  of  the 
PEM  fuel  cell  have  been  developed  by  other  workers 
[1,10,22-24].  The  complexity  of  the  model  varies  consider¬ 
ably.  A  typical  model  proposed  by  Lee  et  al.  [24]  presents 
a  comprehensive  empirical  relation  between  voltage  and 
current,  i.e. 

V  =  Vo  —  b In  i  —  Ri  —  m  exp (ni)  —  bln  ^  (17) 

where  b,  R,  m  and  n  are  empirical  parameters  to  be  deter¬ 
mined  experimentally.  The  first  term  on  the  right-hand  side 
of  Eq.  (17)  is  the  open-circuit  voltage,  usually  1.2  V  for  a  hy¬ 
drogen  fuel  cell.  The  remaining  four  terms  on  the  right-hand 
side  of  Eq.  (17)  account  for  different  kinds  of  polarization 
(or  irreversibilities)  of  the  cell.  The  second  term  is  called  the 
‘activation  polarization’,  caused  mainly  by  the  slowness  of 
the  reaction  on  the  electrode  surface.  The  third  term  is  the 
‘ohmic  polarization’  due  to  the  resistance  to  the  flow  of  ions 
through  the  electrolyte.  These  first  three  terms  are  essen¬ 
tially  the  same  for  the  different  models  proposed  before  and 
can  provide  a  good  prediction  of  the  performance  when  the 
current  density  is  lower  than  the  limiting  value.  The  last  two 
terms,  usually  the  concentration  polarization,  are  proposed 
particularly  for  the  performance  when  the  current  density  is 
close  to  the  limiting  current  density  and  the  voltage  experi¬ 
ences  a  rapid  drop.  Recently,  a  large  effort  has  been  devoted 
to  discussing  the  mathematical  form  of  the  fuel  cell  perfor¬ 
mance  in  this  high  current  density  regime.  There  is,  however, 
a  lack  of  physical  justification  of  the  corresponding  math¬ 
ematical  expression.  For  example,  the  exponential  term  of 
Eq.  (17),  first  proposed  by  Kim  et  al.  [23],  turns  out  to  serve 
merely  as  a  tool  for  curve  fitting  so  that  the  parameters  m 
and  n  have  to  be  recalculated  for  any  change  in  the  operation 
conditions.  The  logarithmic  term  proposed  by  Lee  et  al.  [24] 
is  a  form  of  the  Nernst  equation;  it  is  of  a  more  clear  physical 
sense  and  accounts  for  the  influence  of  the  oxygen  fraction 
on  the  voltage  drop.  Namely,  a  smaller  oxygen  mass  fraction 
leads  to  a  larger  voltage  drop,  and  this  is  supported  by  the 
present  results.  Actually,  the  present  computational  results 
indicate  that  the  level  of  polarization,  or  the  voltage  drop,  in 
the  high  current  density  regime  corresponds  largely  to  the 
transport  of  oxygen  at  the  cathode,  and  this,  as  shown  by  the 
results  discussed  in  the  previous  section,  is  related  closely 
to  the  porosity  of  the  GDL.  A  limiting  case  of  this  scenario 
is  that  as  soon  as  the  oxygen  is  consumed  completely,  the 
voltage  shall  drop  to  zero  because  the  fuel  cell  is  expected 


Fig.  7.  Polarization  and  power  density  curves  corresponding  to  four 
different  porosity  models. 

to  cease  functioning  due  to  the  lack  of  oxygen.  This  occurs 
when  the  system  approaches  the  limiting  current. 

The  above  scenario  is  examined  in  the  present  analysis  by 
considering  the  four  porosity  models  and  the  fuel  cell  perfor¬ 
mance  is  discussed  with  the  aid  of  polarization  curves.  For 
the  base  case  conditions  [20],  there  are  various  relationships 
between  any  two  of  the  following  parameters:  the  external 
load,  the  cell  voltage,  the  cell  operating  current  density,  the 
surface  overpotential  across  the  cathode  catalyst  layer.  When 
the  surface  overpotential  is  taken  as  the  governing  parameter 
of  the  system  [20],  the  polarization  curve  can  be  described 
parametrically  by  the  relation  between  V(rj)  and  where 
V(r])  =  Vo  —  A<P(rj)  —  t)  is  the  cell  voltage  and  Vo  is  the 
open  circuit  voltage.  The  power  density  curves  also  can  be 
described  parametrically  by  the  relation  between  PW(rj)  and 
where  PW(r/)  =  V(rj)  x  i(rj)  is  the  power  density.  The 
curves  of  polarization  and  power  density  corresponding  to 
the  four  porosity  models  are  shown  in  Fig.  7,  in  which  a 
comparison  between  the  polarization  curve  of  the  present 
model  1  and  that  of  the  analytical  work  of  Gurau  et  al.  [20] 
for  a  uniform  porosity  £dW  =  0.4  is  also  shown.  It  is  found 
that  the  present  numerical  results  are  in  excellent  agreement 
with  the  analytical  results  of  Gurau  et  al.  [20],  which  con¬ 
firms  that  the  present  numerical  approach  has  been  correctly 
implemented. 

The  polarization  curves  consist  of  three  regions  of  differ¬ 
ent  voltage  drop,  which  is  typical  for  a  PEM  fuel  cell.  They 
are:  (i)  the  region  on  the  left  (first  sharp  drop)  corresponds 
to  the  activation  polarization,  accounted  for  by  the  second 
term  of  Eq.  (17);  (ii)  the  region  in  the  middle  (decreasing 
linearly)  corresponds  to  the  ohmic  polarization,  accounted 
for  by  the  third  term  of  Eq.  (17);  and  (iii)  the  region  on  the 
right  (second  sharp  drop)  corresponds  to  the  concentration 
polarization,  accounted  for  by  the  last  two  terms  of  Eq.  (17). 
Note  that,  the  polarization  curves  of  the  four  porosity  mod- 
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els  overlap  in  the  region  of  low-to-medium  current  density 
(the  linear  regime),  while  at  higher  current  density  close  to 
the  limiting  value,  differences  between  the  four  models  be¬ 
come  significant.  This  implies  that  the  porosity  of  the  GDL 
exerts  no  influence  on  fuel  cell  performance  when  the  cur¬ 
rent  density  is  of  a  medium  or  low  value,  but  has  a  signifi¬ 
cant  when  the  current  density  is  close  to  the  limiting  value. 
The  explanation  suggested  by  the  present  results  is  that  a 
larger  porosity  in  the  GDL  results  in  a  larger  oxygen  trans¬ 
fer  from  the  gas  channel  to  the  catalyst  layer,  and  thus  a 
larger  current  density  is  generated.  Thus,  for  example,  for 
the  model  of  the  highest  porosity,  i.e.  model  3,  the  limiting 
current  density  is  the  largest.  This  finding  implies  that  the 
polarization  curve  in  the  regime  close  to  the  limiting  current 
density  is  governed  by  the  performance  of  the  mass  transfer 
across  both  the  GDL  and  the  catalyst  layer,  as  commonly 
proposed  in  previous  studies  [23-25].  Similar  phenomena 
were  also  found  by  Rho  et  al.  [26],  who  indicated  that  the 
drastic  drop  of  voltage  is  due  to  the  transport  limitation  of 
oxygen  through  the  GDL,  and  a  mixed  gas  with  a  higher 
oxygen  leads  to  a  higher  limiting  current  density.  Rho  et  al. 
[27]  further  pointed  out  that  the  departure  of  cell  potential 
from  the  linear  regime  is  caused  by  partial  flooding  of  the 
GDL,  which  decreases  the  oxygen  diffusion  in  the  GDL  and 
thus  the  voltage  of  the  cell. 

It  is  also  seen  from  Fig.  7  that  the  power  density  increases 
virtually  linearly  with  current  density.  When  the  latter  ap¬ 
proaches  the  limiting  value,  the  power  density  reaches  a 
maximum  and  then  decreases  sharply  because  of  the  deple¬ 
tion  of  oxygen  in  the  catalyst  layer.  The  maximum  power 
density  increases  and  corresponds  to  a  higher  current  den¬ 
sity  as  the  porosity  of  the  GDL  is  increased.  Consequently, 
a  fuel  cell  can  generally  deliver  better  performance  when  a 
GDL  of  larger  porosity  is  used. 

5.  Concluding  remarks 

The  influence  of  the  porosity  of  the  GDL  on  the  perfor¬ 
mance  of  a  PEM  fuel  cell  has  been  investigated.  Four  poros¬ 
ity  models  which  represents  different  continuous  functions 
of  position  are  considered  so  that  the  porosity  influence  can 
be  studied  in  a  systematic  way.  Due  to  the  complexity  caused 
by  the  porosity  models  considered,  the  resultant  simultane¬ 
ous  differential  equations  are  solved  numerically.  Results 
corresponding  to  the  uniform  porosity  model  are  in  excel¬ 
lent  agreement  with  the  analytical  results  of  Gurau  et  al. 
[20],  and  thus  confirm  the  accuracy  of  the  present  analyses. 
It  is  also  found  that  the  influence  of  the  porosity  change  in 
the  GDL  generally  reflects  the  oxygen  transfer  from  the  gas 
channel  to  the  GDL  and  catalyst  layer,  i.e.  greater  porosity 
leads  to  higher  consumption  of  oxygen  such  that  higher  cur¬ 
rent  density  is  generated  and  a  better  fuel  cell  performance 
is  obtained.  The  present  analyses  also  serve  to  explain  the 
physical  meaning  of  the  performance  model  in  the  high  cur¬ 
rent  density  regime,  which  has  been  an  issue  discussed  in 


many  previous  publications.  In  summary,  the  following  con¬ 
clusions  have  been  drawn. 

1.  A  higher  averaged  porosity  of  the  GDL  results  in  a 
smaller  oxygen  mass  fraction  in  the  gas  channel  because 
a  higher  consumption  of  oxygen  in  the  catalyst  layer 
occurs.  Consequently,  a  fuel  cell  with  a  higher  surface 
overpotential,  r),  may  require  a  GDL  of  larger  porosity 
so  that  a  larger  quantity  of  oxygen  can  be  provided  from 
the  gas  channel  to  the  catalyst.  This  is  because,  due  to 
the  higher  more  oxygen  is  required  to  support  the 
active  electrochemical  reaction. 

2.  A  continuous  change  in  the  porosity  of  the  GDL  may  be 
necessary  to  accommodate  the  fact  that  the  oxygen  mass 
fraction  is  not  uniform  across  the  GDL,  especially  in  the 
region  adjacent  to  the  interface  between  the  GDL  and 
the  catalyst  layer.  This  effect  is  less  significant  when  r\  is 
large  because  the  electrochemical  reaction  is  then  more 
active  and  more  oxygen  is  consumed,  which  leads  to  a 
smaller  difference  oxygen  in  the  oxygen  mass  fraction 
across  the  GDL. 

3.  Across  the  catalyst  layer,  the  current  density  increases 
rapidly  from  zero  at  the  interface  between  GDL  and  cata¬ 
lyst  to  the  limiting  value  at  the  interface  between  catalyst 
and  membrane.  The  limiting  current  density  increases 
with  r]  because  a  higher  r]  leads  to  a  more  active  electro¬ 
chemical  reaction  as  well  as  a  higher  current.  The  current 
density  also  increases  with  the  averaged  porosity  of  the 
GDL  because  of  the  larger  consumption  of  oxygen  in  the 
catalyst  layer. 

4.  A  higher  porosity  in  the  GDL  leads  to  a  larger  current 
density  and  therefore  a  higher  loss  of  the  membrane  phase 
potential,  which  also  corresponds  to  a  larger  overpotential 
and  a  better  performance  of  the  fuel  cell. 

5.  A  change  in  the  porosity  of  the  GDL  has  virtually  no  in¬ 
fluence  on  the  polarization  level  when  the  current  density 
is  at  a  medium  or  low  value.  On  the  other  hand,  it  has 
a  significant  effect  on  the  polarization  when  the  current 
density  is  close  to  the  limiting  value.  This  is  obviously 
due  to  the  fact  that  a  larger  porosity  in  the  GDL  results 
in  a  larger  oxygen  transfer  from  the  gas  channel  to  the 
catalyst  layer,  and  thus  a  larger  current  density  is  gen¬ 
erated.  This  result  supports  the  fact  that  the  polarization 
curve  in  the  regime  close  to  the  limiting  current  den¬ 
sity  is  governed  by  the  performance  of  the  mass  transfer 
across  both  the  GDL  and  the  catalyst  layer,  as  proposed 
in  previous  studies  [22-24] . 

Although  the  present  results  indicate  coherently  that  a 
larger  porosity  of  the  GDL  leads  to  a  higher  current  den¬ 
sity,  and  to  some  extent  a  better  performance  of  the  fuel 
cell,  a  high  porosity  is  accompanied  by  water  flooding  in  the 
GDL  and  this  decreases  markedly  the  cell  voltage  [16].  Be¬ 
sides,  as  concluded  by  the  present  paper,  a  fuel  cell  with  a 
higher  surface  overpotential,  r\,  may  require  a  GDL  of  greater 
porosity  because  more  oxygen  is  necessary  to  support  the 
active  electrochemical  reaction.  Gurau  et  al.  [20]  indicated, 
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however,  that  a  higher  rj  induces  a  higher  ohmic  resistance 
in  the  catalyst  layer,  which  in  turn  decreases  the  fuel  cell 
voltage  under  moderate  current  density  because  ohmic  re¬ 
sistance  dominates  the  cell  performance. 

It  is  also  noted  that  optimization  of  the  GDL  depends  not 
only  on  porosity  morphology,  but  also  on  the  thickness  and 
the  content  of  the  GDL.  Giorgi  et  al.  [18]  have  shown  that  a 
change  in  the  content  of  PTFE  of  GDL  results  in  a  change 
in  porosity,  and  thus  affects  the  cell  performance,  although 
the  influence  was  insignificant  in  their  experiments.  They 
further  pointed  out  that  the  thickness  of  the  GDL  has  a  more 
substantial  effect,  i.e.  a  thick  GDL  may  decrease  the  cell 
performance  because  the  oxygen  path  is  too  long;  a  thin 
GDL  may  also  decrease  the  cell  performance  because  of  the 
contact  resistance  at  the  GDL  I  catalyst  interface. 

In  summary,  there  is  an  optimum  porosity  level  for  the 
GDL  so  that  the  fuel  cell  has  the  best  performance.  This 
optimization  is  beyond  the  scope  of  present  study,  which 
can  only  be  implemented  when  the  whole  fuel  cell  system 
is  considered  and  a  more  thorough  mathematical  model  is 
analyzed. 
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